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ABSTRACT 


The nonlinear evolution of unstable sound waves in a uranium 
plasma has been calculated using a multiple time-scale asymptotic 
expansion scheme. The fluid equations used include the fission power 
density, radiation diffusion, and the effects of the changing degree of 
ionization of the uranium atoms. The nonlinear growth of unstable waves 
is shown to be limited by mode coupling to shorter wavelength waves which 
are damped by radiation diffusion. This mechanism limits the wave pres- 
sure fluctuations to values of order 6P/P ~ 10“ 5 .in the plasma of a 
typical gas-core nuclear rocket engine. The instability is thus not 
expected to present a control problem for this engine. 
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1. INTRODUCTION 


The objective of this work is to analyze the nature of instabilities 

in fission reactors in which the uranium operates at a temperature sufficiently 
high so as to become ionized. The principal impetus for study of such uranium 
plasma reactors has derived from the concept of a high specific impulse nuclear 
gas-cored rocket engine.' ■ * However, fissioning plasma reactors may also lead 
to development of a new class of magnetohydrodynamic generators' 3 *, and perhaps 
may even have application to the direct nuclear excitation of lasers. 


The equations describing a uranium plasma involve the coupling of the 

to the neutron transport and also to the radiation trans 
£ort in the system. In such a reactor, the plasma would be enclosed (or 
nearly enclosed) by a neutron moderator - whi ch is important in that it returns 
the thermal i zed neutrons back to the plasma core - while radiation transport 
dominates the transport of energy within the plasma. As a theoretical structure 
en, the coupled systen of equations can be charted symbolically as: 



Figure 1.1 


The specific set of equations involved are listed as equations (2.24)-(2.28) 
the next section of this report. 


in 


The major application for the fissioning uranium plasma reactor (on 

w ich significant effort has been expended) is that of the gas-cored nuclear 
rocket engine shown below. 






2 . 




Figure 1.2 


In general, this system consists of a chamber in which a fissioning uranium 
plasma (at ~ 5 * 10 4 °K) is maintained in the central region and radiatively 
heats the propellant hydrogen gas which flows around the main core and is 
subsequently ejected from the rocket nozzle. The chamber is surrounded by a 
moderator which slows down the fission neutrons and returns them to the 
uranium core. 

It is clear that if instabilities occur in this device, leading to 
either a fine-scale plasma turbulence or larger scale oscillations, they could 
well have undesirable results — such as enhancing the uranium fuel ejection 
from the rocket. An understanding of these processes (and their possible 
controllability) is of basic importance for the feasibility of constructing 
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such plasma-cored rocket engines. 

In a more general context, it should also be noted that the direct 
excitation of plasma motion and plasma waves by a fissioning plasma reactor 
may have useful applications to other energy conversion devices such as MHD 
generators . 

One source of instability is the growth of sound waves^ 4,5 ^ in the 
plasma due to the fission power density, P^. . The instability mechanism 

involved can be understood as follows: Consider a standing sound wave in a 

bounded region of fissioning plasma with a constant background thermal neutron 
density. In the wave compressions the fission power density, P f . , increases 
due to the increased uranium density, and in the rarefactions P-. decreases. 
This results in an increased pressure gradient associated with the wave which 
in turn leads to a transfer of fission power to the wave. It occurs because 
the wave compression tends to expand more rapidly than it was compressed. How- 
ever, competing with this is the fact that radiation tends to transport the 
extra thermal energy out of the wave compressions. Radiation diffusion smooths 
out the temperature fluctuations of waves more rapidly the shorter their wave- 
length. This results in a critical wavelength, A . t » below which waves are 
stable and above which they are unstable. 

We have calculated the nonlinear evolution of such unstable waves 
to determine their limiting amplitude. The mechanism that finally limits 
their amplitude is the mode coupling of unstable waves in the wavelength range 
A > A cn - t to radiation damped waves in the range A < A cr j t 

Our conclusion is that this mechanism "turns off" the instability 
after the wave amplitudes have reached relatively small values for most sys- 
tems of interest. Relative pressure fluctuations of order 6P/P ~ 10 -7 to 
10 -5 are calculated to occur for a slab of uranium plasma with parameters cor- 
responding to the NASA-Lewis design engine. Consequently this instability is 
not expected to present a control problem for this engine. The instability 
would only become a problem for devices of a size much larger than those cur- 
rently contemplated. 



A second mechanism for instability derives from the fact that the 
fission power density in the reactor cavity depends on the relative fluid 
velocity between the neutron gas and the uranium plasma. For low neutron 
thermal energies (below 1 eV which is the lowest neutron energy for resonant 
capture by U235 ) the neutron induced fission power increases as the relative 
fluid velocity increases. This "Doppler effect" enhances the fission power 
in some parts of a sound wave and not others and can give rise to growth of 
the wave. This process is analysed in Chapter 7. Our conclusion is that 
nonlinear effects also limit the amplitude of waves in this case to small 
values which are not expected to cause problems for the gas-core engine. 
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5 . 


2. DERIVATION OF FLUID EQUATIONS FOR FISSIONING U-PLASMA 


In this section the fluid equations for a uranium plasma will he 
derived from the equations for the particle distribution functions. This has 
the advantage that one can see physically the assumptions that are involved 
in the fluid description and where future generalizations might be made. 

Define a distribution function F so that 

F(a, v_, x_, t) d_v da = (average number of U atoms or ions (2.1) 

in dv_ in states a, da. ) 

If the states labelled a are discrete for part of the a range, then F 
involves a series of 6-functions in a. The corresponding velocity distri- 
bution function for U atoms in any charge state is 

f = / F da , (2.2) 

and the total number density for all U atoms and ions is 

N u = / f d v = 2 , (2.3) 

where i goes over all the partial densities for ions that have lost i 
electrons . 


Similarly the fluid velocity V u and pressure tensor P U _I_ 
for the U-component of the plasma are written 

(2.4) 

(2.5) 

where we have assumed a local velocity distribution close to a Maxwellian 


M, - /“vvf . 

P I ’ V T I * M / dV f (V - V,,) (V - Vy) , 
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with a temperature T. Viscosity, which derives from the off-diagonal parts 
of P u , will be neglected. 

Now the transport equation is 


3F 

3t 


+ 


V 


3F 

3x 


( 2 . 6 ) 


where C u dx_ dv_ da = rate of change of U-atoms in dx_ d_v da due to atomic 
processes, fission processes, and interaction with the radiation field (i.e . , 
C includes all microscopic effects). 


Taking moments J d_v da and J _v d_v da of (2.6) gives the mass and 
momentum transport equations. 


3N 


u * lr ' < N uV ■ f c u d * d “ * 


3t 


(2.7) 


It < N U V * lx <V T > + 1? ♦ < N AV ■ / c u i d x d “ 


(2.8) 


For the energy transport equation we operate on (2.6) with 


J dv_ da ^ M (v_ - v^) 2 + e(a) .... 


'thermal energy 
density excluding 
.translational K.E.) 


( “internal energy" of ^ 
ionization excitation 
and radiation 7 


Wri ti ng 





It 


(2.9) 



the resulting energy transport equation becomes. 


d ( 3 


at It V T ) + It (N u e) + Ir • + 1 V T I? • \ + la ' ^ 


= Jdv da C u \ M (v - Vj 2 + e 


( 2 . 10 ) 


where 


d_ = 3_ . » 9_ 

dt . 3t — u 3x 


( 2 . 11 ) 


In order to obtain a closed set of 1-fluid equations we must next 
consider the meaning of the collision terms in (2.7) - (2.10) and the role 
of the electrons. 


Number Density Equation 

There is no local source of U atoms, and we shall neglect the local 
sink due to fission. We also neglect the accumulation of a spectrum of fis- 
sion fragments which would form a new ion species. Under these circumstances 
J dv_ da = 0, so that equation (2.7) gives 


3N 

3t 


u 


+ 


3 

3x 


ayy - o 


( 2 . 12 ) 


for the continuity equation for U number density. 


However, for the electrons there are sources and sinks due to the 


rapid atomic processes of recombination and ionization. These atomic processes 
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occur on time scales much less than those associated with the fluid dynamics. 
Thus a value for the electron number density N g will be assumed corresponding 
to local thermodynamic equilibrium 

N e = N U Z (2.13) 

where Z is the mean charge state of the U atoms and must be calculated from 
the Saha formulas. 

Momentum Equation 


Fission does not change the average local momentum density of the 

plasma. The electrons and ions however couple strongly through collisions 

so that V = V, T = T = T, and Z m f v dv C = 0. 

-e -ire u „ • J — — 

e,i J 

Thus adding (2.8) to it's counterpart for the electrons and treating 
m/M u << 1 gives 


M 


3t 


< N uV 


9 

9x 


(N, 


N e )KT 


+ M 


9x 


(NJWi) = 0 


(2.14) 


Note, we could include a photon pressure in (2.14) but have not done so 
at this point. 


Energy Equation 

Now the local energy density (excluding only translational kinetic 

energy) is 


/< 


dv da 


(e 


exct. 


+ e 


lomz 


e Rad^ 


M(v - 


V ) 2 


(2.15) 


where the quantity e in (2.10) has been split up into its contributions from 
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internal energy of excitation, ionization and radiation. The ionization energy 
is viewed as similar to internal excitatio-n energy for the atoms. The various 
atomic processes involved undergo energy interchanges with the thermal velocity 
component and not the ordered kinetic energy density of the plasma, which is why 
we chose to take the above moment for the energy equation. The quantity (2.15) 
changes due to collisions, but when summed over atoms and electrons its remain- 
ing changes (due to C) are those due to fission and radiation, i.e. 3 


l 

e,u 


d_v 


dv da C 


'exct 


+ £. 


horn z 


e Rad^ 


M( v - 


Xu ) 2 


= P 


fi ss 


+ P 


v 


(2.16) 


where 


P fiss = fission Power density in frame of the U-plasma. 

P v = nett local radiation energy absorption less emission 


(2.17) 


Adding the electron version of (2.10) to the ion energy equation 
(2.10) then gives the total energy transport equation for the fluid. 


d_ 

dt 


i V 1 * 1 + z >] + W (N u^ • (i v> + 1 o + z) 




A, 


+ 


3 

3x 


* 3e> = p fiss + p v 


(2.18) 


The thermal conduction terms and will be neglected since the most 
important energy transport is due to radiation as expressed by P The 



quantities and are next required in order to complete the fluid 

equati ons . 

Fission Power Density, P-. 

J ’ fi ss 


Suppose v Q is the average thermal velocity of the neutrons, W Q 
their number density, and a the cross-section for neutron induced fission 
events in which an energy Q is released in fission fragments. The fission 
power density then follows from 

I total path traveled \ n 
(W 0 v 0 dx_ dt Q) y by n Q 1 s in d><_ dt J ^ 

(1/N a) (m.f.p. for fission) 

i.e.j 


p fiss dt = 


P fiss 


% N u v o a 0 


(2.19) 


We have assumed c to be independent of the relative velocity 
between the neutrons and U nuclei and taken a value for it averaged over 
the neutron distribution function, i.e., a = a(v Q ,T). A more precise form 
for Pf.} ss would be 

P fiss = // d *C d ^ 1*0 - f o<V - *ol> (2 - 20) 

where Q 1 is the fission kinetic energy release in the rest frame of the 
plasma, f and f y are velocity distribution functions for the neutrons and 
U nuclei. The simpler model (2.19) will be used for most purposes. 



Radiation Transport term, P v 


In the interior region of the U plasma the radiation absorption 
length is usually much less than the gradient scales of the plasma, i . c . , 
the medium is optically thick. Under these circumstances the absorption 
coefficient can be approximated by the Rosseland mean value, K R , given by 


16 oT 3 
R " 3 k R 


(2.21) 


where a is the Stefan-Bol tzmann constant, T the temperature, and k R the 
Rosseland mean opacity. The radiation energy flux is then given by 


S = - K r V_ T (2.22) 

and the divergence of this, - V_* S, gives the local radiation power, P^, 

l.o.. , 


P v = V_ • [ K r v_ T] . (2.23) 

The quantity K R is thus an effective thermal conductivity due to radiation 
transport. 


Summary of Fluid Equations 

Equations (2.12), (2.13), (2.18), (2.19) and (2.23) now form a 
complete set of fluid equations. Since they are no longer necessary, we shall 
drop the subscripts u on the total ion number density N u and fluid velocity 
V , and also drop the bars over F and T. The resulting system is, 


3N_ 8_ 

at ax 


• (NV) = 0 


(2.24) 
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a -i pi 

Bt" + — * 3x" — + MN" W 


NKT (1 + Z) 


(2.25) 


d_ 

dt 


| NKT (1 + Z) 


ft (Ne) + • (NeV) + | (1 + Z) NKT 


3 

3x 


(2.26) 


P fiss +V -' (Kr^ T) 


flss 


= N v n 0 Q 
oo ^ 


w 1 u 

I'n n 


16 aT 3 
3 ''R 


(2.27) 

(2.28) 


where Z, e, must be calculated from the Saha relations for a thermal 
equilibrium U plasma. 

Comments on the Energy Equation and Specific Heats 

If we write the total "internal energy" (which includes thermal 
kinetic energy) as 

E = | NKT (1 + Z) + Ne (2.29) 

then equation (2.26) can also be written 


f +E | - V + = 0 , (2.30) 

where the pressure P = (1 + Z)NKT. This is the usual form of the energy 
equation (see Chapman & Cowling, Mathematical Theory of Nonuniform Gases, 
p.52, Eq. 4). 
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One can also derive an equation that expresses local energy con- 
servation for the total energy. Multiplying (2.25) by MNV, adding the 
result to (2.26) and using (2.24) leads to 


~ \ MNV 2 + | NKT (1 + Z) + Ne 


3x 


j MNV 2 + | V NKT (1 + Z) + NVe - K 


3T 
R 3x 


= P 


fi ss 


(2.31) 


This shows that the total energy density in the medium is 


‘tot. 


MNV 2 + E 


(2.32) 


The specific heats can be expressed in terms of these quantities 
(see Landau and Lifshitz, Stat. Mech. p. 45, or Chapman and Cowling, p.39). 
Thus 


3E\ r _ / aCE + Pv] 

3T J ’ L p " \ 3T 


where v is the gas volume. From (2.29) we see 


C v = NK 


L f (1 + z + TZ t ) + ^ 


(2.33) 


(2.34) 


The effective number of degrees of freedom, n, in the plasma can 
be defined by 


E 


j n NKT , 


(2.35) 




so from 


However, 
quanti ty 


I 



(2.29) 

n = 3 (1 + Z) + —■ . 

since n is a function of T, this is not a particularly useful 
in that the quantity y s C p /C v ^ (1 + ^-). 

The first term, dE/dt, of (2.30) can also be written C y dT/dt. 


14 . 


(2.36) 



3. NUMERICAL RESULTS FOR THE AVERAGE IONIZATION ENERGY e AND AVERAGE 
CHARGE" NUMBER Z OF THE U- IONS 


The functions e and Z in this section are derived from some com- 
puter results kindly supplied by Dr. R. T. Schneider^. They are based on 
the UPLAZ-2 program (developed by the group in the Department of Nuclear 
Engineering Sciences, University of Florida) which, calculates the composition 
of uranium plasma as a function of plasma pressure and temperature assuming 
local thermodynamic equilibrium. 


The average ionization energy e is defined as 


e 


Z N 


( 1 ) 


Z N 


ITT 

u 


(3.1) 


where the sum goes over all ion charge states, and e- is the total energy 
required to produce an ion which has lost i electrons (including energy 
shifts due to the density effect), and is the number of ions in 
ionization state i per cm 3 . 

The energy e. can be written as a sum due to successive ioniza- 
tions, 

£ i = ( E 1 " AE i) + ( E 2 " AE 2^ + •••( E i " AE i)> £ 0 = °» (3.2) 

where E-j is the energy required for the process 

•> U + + e , 

and E 2 for U -»• U + e, etc., and the AE's give the decrease in these 
energies due to the density effect. 
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The average charge number Z of the ions is also required and is 

defined by 


Z 


Z Zl 


,(i) 


E N 


(i) 


(3.3) 


The quantities N^ 1 ^ , AE i , f ZZ- N^ 1 ^ and (Z + 1)N U = N t otal’ 
were supplied as functions of T and ^ (total pressure) by UPLAZ-2, 
and after some computation the tables and graphs 3.1 and 3.2 of 
e(T, P tQt ) and Z(T, P tot ) were obtained. 




P. rt . = 10 Atm. 
tot 

P, . = 100 Atm. 
tot 

P, ,= 500 Atm. 
tot 

P. , = 1000 Atm. 
tot 

T 

Z 

Z 

Z 

Z 

120,000 

9.552 

8.273 

7.053 

6.850 

100,000 

7.951 

6.490 

6.058 

6.058 

80,000 

6.220 

5.843 

5.375 

4.993 


5.421 

4.810 

4.231 

4.139 

40,000 

3.692 

3.463 

3.168 

3.159 

20,000 

1.913 

1.368 

.848 

.682 

8,000 

.159 

.057 

. 

.031 

.025 

_ . ... 

5,000 

.006 

.002 

.001 

.001 


T 

P. .= 10 Atm. 
tot 

e eV 

P. 100 Atm. 

tot 

£ 

P, ,= 500 Atm. 
tot 

£ 

P. 1000 Atm. 

uOt 

£ 

120,000 

552 

383 

261 

221 

100,000 

363 

210 

167 

154 

80,000 

188 

153 

117 

93 

60,000 

129 

89.4 

63.5 

55 

40,000 

55 

42.9 

31.2 

28.2 

20,000 

16.6 

9.1 

4.6 

3.5 

8,000 

.92 

.3! 

.17 

1 .3 

5,000 

4 • 1 0 - 2 

10" 2 

7* 10" 3 

5 • 10” 3 
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4. DISPERSION RELATIONS FOR SOUND WAVES IN FISSIONING URANIUM PLASMA 


The first step in our analysis is to derive the dispersion relations 
for sound waves in a fissioning U-plasma by linearizing the fluid equations in 
the wave amplitudes. The plasma in which the waves propagate is assumed to 
be stationary in time and approximately spatially homogeneous on the wavelength 
scale of the waves being considered. 

The fluid equations are (dropping subscripts u and bars on e and 

Z), 


3N 

3t 


+ 


3 

3x 


• (NV) = 0 


(4.1) 


3V 

at 


+ V-| E V + 1 F | S [wcr(ltZ)] - o 


(4.2) 


i_ + v . i- 

at - 3x 


j NKT ( 1 + Z) + Ne 


+ | (1 + Z) NKT |^ • V 


+ V * P fiss + - ’ (K R^ T) 


(4.3) 


where 


f P-. = N ALv cr Q 

fiss 00 s 


Z, c, are given fns of T, N 


(4.4) 


Now write the fluid variables as 



N = 


<o + N l 


V = 0 + V 1 


T = T o +T l 


Z = Z o + N 1 Z N + T 1 Z T 


£ = £ o + N 1 £ N + T 1 £ T 


(4.5) 


x 

where = (3Z/3N)^_^ y_-|. etc., and the equilibrium state is denoted by 

subscripts o and perturbati 8ns by subscripts 1. 


The linearized equations then follow as. 


3N 1 3 

dT + H ok * 


= 0 


(4.6) 


3V, 1 3 

3t + MN^ 3x 


N l KT o (1 + Z o + W + N o KT l + Z 0 + T o Z T>J = 0 < 4 - 7 > 


3_ 3_ 

2 at 


’n,KT (1 + 1 + N Z w ) + N KT, (1 + 1 + T Z T ) 
L I o o on o 1 o o T 


+ 


N n e„ + N„N-,e M + N 0 T 1 £ t ] + N_e_ ^rr • V- 


n -i C. ' H.l'll 

at l i o o i 


'o^o ax -i 


+ I 0 + Z o> N o KT o It • = 


N 


1 


N q P fiss,0 + K R0 V2 T 1 


(4.8) 


Next, consider Fourier modes 


(N,, V_-| , T -j ) -»■ (N-j , V,, T-j ) 
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Thus 


wN-, + N 0 J< ' i-j = 0 


1 


- ^ + k ( AN-j + BT 1 ) = 0 

o 


It - 0)( r A + e 0 + N 0 £ n ) + 


iP., 


fiss,o 


+ T ] - to(| B + N 0 e T ) - i k 2 K 


V. L 


I 0 + z o) n 0 kt 0 + n 0 £ 0 


RO 


k • V , = 0 


(4.9) 


(4.10) 


(4.11) 


where 


A - KT 0 (1 ♦ Z 0 + N 0 Z n ) 


B = N q K (1 ♦ Z 0 + T 0 Z t ) 


(4.12) 


Now from (4.9) 


!•> = — k • V, 

1 to — — 1 


so (4.10 and (4.11), become. 


B T n 


■ w -i + Mw - -i ) + - wr~ 0 


(4.13) 


k-V, j - N o (| A + e 0 + N qEn ) + + f (1 + Z 0 )N 0 KT 0 + 

+ T 1 | - “(f B + V T ) - ik 2 K R0 j = 0 . 


(4.14) 


Eliminate T-j from (4.13) 


- to V ] 



(J<“ V-] ) 


k JL 

- m. 


(k * ^ ) rT - 


= 0 


(4.15) 


The modes all have along i£ so that we can take k_* and cancel _k*_V-j 
thereby obtaining a dispersion relation 


- co + 


Ak : 

Mco 


) {- B + N o £ l) ‘ lk2K R0 } 2 " KN^~ {“ N o (I A + £ o + N o £ n) 


iP.e.t 


+ | (1 + Z )N KT + 

to 2 v 0 ' 0 0 


I £ \ 
0 0^ 


= 0 


(4.16) 


This dispersion relation is a cubic in co and can be written, 


to' 


^ B + N Q £y) + to 2 ik 2 K ro + co I - 


Ak 2 

M" Vt 


■ k 2 B / , M \ k 2 B 

M ( £ o N o £ N ) “ MN 


0 


| (1 + ZJN KT + N^c. 
2 0 0 0 0 0 


ll 


(4.17) 


Aik 1 * y , 2n lP fiss,0 l _ n 

IT K R0 + k B MN " 0 • 


Free Sound Wave Case 


If K rq = 0 and P fiss Q = 0, (4.17) becomes the .relation for a 
simple sound wave (except for the peculiarity of a changing degree of ioni- 
zation in the U-plasma as the temperature fluctuates). In this limit we 
find 


= ± k V $ , 


co 


(4.18) 




A |>! 

M N o £ T 


M ( £ o + N o £ N ) + M 


! KT o + Z o> +£ o 


]}* 


|I b + Vt} % 


where V is the sound speed. Note also that dw/dk = co/k in this case. 


If we use A and B given by (4.12), the sound speed can be 
written in the form 


V s = 


5KT 

0 
3M 


3 /n . , x 3 N o £ N , 3e T ^ + Z o + N oVi 
2 U V ~ 5KT Q + 5K(1 + Z Q + T^F j 


h 


2 K(1 + Z 0 + T q Z t ) 


General Case including Fission and Radiation 

Two fairly simple limits can be extracted from the cubic (4.17). 
The first corresponds to the case in which the competing effects of radi- 
ation and fission nearly balance. The result is a growth or damping of the 
wave that is small. Thus we can write 

to = oo n + oo T = kV_ + W T 
R I _ s I 


in (4.17) and treat toj as small (this gives an exact condition for the 
critical wavenumber for which ojj changes sign). The cubic can be rewritten 



Vt 


(jJ (oo 2 


k 2 V 2 ) + w 2 



Aik v K R0 

M 


k 2 BiP 


MN 


f i s s , o _ Q 


24 


(4.19) 


(4.20) 


(4.21) 


(4.22) 
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Linearizing in ojj gives 
^BP 


co 


I 


■f i SS , 0 i, 2w2i/ . 

”” k V s K R0 

0 

2N 

3 V s< B + -T- £ I > 


Ak 2 K 


RO 


(4.23) 


The system is unstable if the numerator of (4.23) is positive. 
There are two cases for which this is theoretically possible, namely 


r 


k| < k 


cri t 


BP,, 
fi ss ,o 


h 


ho m o [% 


for V* > i 


k >0 for V 2 < A/M 


(4.24) 


The inequality V 2 > A/M can be written 


N e N o 

1 + z o> -kt^IVn 


(4.25) 


The first range of unstable waves with wavelengths greater than a critical 
wavelength is similar to that found by Becker and McNeil. The expression 
(4.24) for the critical wavenumber, k at which waves go from stable to 
unstable is exact, but the growth rate coj for marginally unstable waves is an 
approximation which is valid for | oo-j- 1 « |kV |. It will turn out that 
V 2 > A/M, -c.e., the inequality (4.25) is satisfied, so that only the first 
part of (4.24) is of interest. 
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The second limit that is easily extracted from (4.17) is the 
small k result. We see from (4.17) that to -> 0 as k 0. If we assume for 
small k that to ~ k n (n >0), the various terms in (4.17) are of order k 3n , 
k 2+n , k 2n+2 , k 4 , k 2 . As k 0 only the first and last terms survive. 

Retaining only these two terms and noting i 1 ' 3 = + i)/2, we find the 

small k limit. 


* D k 2 / 3 


2B P 


fi ss ,o 


2N 

3MN(B + e T ) j 


1/3 


(4.26) 


In this case we note dto^/dk = (2/3) to^/k ■+ ® as k 0 . 

The results (4.23), (4.24), and (4.26) can be expressed in terms 
of our original variables by eliminating A and B using (4.12). They become, 


OJt = 


1 j (1 + V Vl> f MN^’ 0 

k2K R0 

N o K 

KT 1 

|/t - M + Z o + N oV_ 

I 

3V 2 

s 

1 + Z o + T I 

3e t ) 

o Z T + 3K j 


(4.27) 


k crit " k " K R0 


KT 

1 - — (1 + Z + NZJ 
MV 2 0 0 N 


3 V 


1 + Z o + T o Z T + K 


T 


(marginal instability) 


2 ^fiss,o 


1/3 


to = 


OJr + lojj 


( /3 + 1 ) k 2/3 


3e 


1 + 


T 


2K(1 + Z 0 T T 0 Z t ) 


(4.28) 


(small k) 
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K (1 + Z 0 + T 0 Z t ) 


p 'i 

fl SS ,0 

M 


^crit ^ 



r KT n 1 


k ro 

[ v s - M + Z o + N oV 1 



(exact) 


(4.29) 


The behavior of |ojj | is shown below in Fig. 4.1. The maximum growth 
rate is approximately the result obtained by setting k = 0 in (4.27). 



Figure 4.1 . Schematic plot of the growth rate of unstable 
sound waves in a fissioning U-plasma as a 
function of wavenumber. 
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5. NONLINEAR BEHAVIOR OF UNSTABLE WAVES AND THEIR LIMITING AMPLITUDE 
IN A U-PLASMA SLAB 

Linear analysis shows that a range of waves with wavenumbers 

k < k . , are unstable. The next question to answer is the limiting ampli- 
cn t 

tude reached as a result of such instability. The principal damping mechanism 
that finally limits the growth of the waves is their mode coupling to that 
part of the spectrum that is radiation damped. This mode coupling should pro- 
duce a steady spectrum of oscillations in a U-plasma, and in this section we 
shall compute this for a slab of plasma of thickness L occupying the space 
0 < x < L. 



Figure 5.1 . 


In order to derive equations for the nonlinear wave amplitudes 
the multiple-time-scale form of perturbation theory^ will be used 
(which is similar to the Bogoliubov-Krylov-Mitropolsky method) through 
terms of second order in wave amplitudes. Thus we express the fluid 
variables as 
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N 

V 

T 

Z 



N + CN, 
o 1 


0 + 

T o + CT 1 

Z o + C ^1 Z N + W 
£q + C(N i£n + T^) 


(5.1) 


and attempt to find an asymptotic solution correct to 0(C 2 ). The parameter 
£ is introduced for bookkeeping purposes and is set equal to unity later. 
The problem is treated as one-dimensional. Second-order derivatives 
Z NN’ Z NT’ etG ' W1 *n be neglected. 


Zeroth-Order System 

The steady-state slab of U-plasma must satisfy the zeroth-order 

equations 



Figure 5.2 
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N KT (1 + Zj = 
oo o 

and 

^fiss ,o + 8x~ ( K R0 

The first complication is now apparent, namely that the zeroth-order slab 
is necessarily inhomogeneous if P-. * 0. The effect of inhomogeneity 

T I jb )U 

will be neglected in this analysis. It will be treated at a later date. 
Equations to 0( C 2 ) 

Substituting (1) into the exact fluid equations and retaining 
terms of 0(c 2 ) leads to a system that can be written as. 


const. = P Q , 


!L 

8x 


= 0 


(5.2) 


(5.3) 


3N, 3V, 

3t~ + N o = ^ 


(5.4) 



3N, 

a 83T +b 





(5.5) 


c 


3N 1 + ri 3T 1 + n 3V 1 

at d at 9 ax 


ca + c 



+ K, 


RO 



(5.6) 


where we have chosen to treat the combined effect of fission and 
radiation on the wave as 0(-£ ) and have neglected the inhomogeneity of 
N , T q , etc-., and defined 
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KT 



a = 


c = 


d = 


(1 + Z n + N a ZJ, b = £ (1 + Z 0 + T 0 z T ) 


'o o 


12” KT 0 <! + Z o + N o Z N> + £ o + N o e N 


i»K |H2 t T„Z t ) + N £ t 
2 o o o T o T 


e " ^fiss.o^o 


g = - (1 + Z ) N KT + N e 
3 2 ' o' o o oo 


(5.7) 


The quadratic terms \p, <J>, A on the right side of (5.4)-(5.6) are given 
by 


* - -lx< N iV 


(5.8) 


4 > = 


MN Q 3x 


N 1 KT 1 0 + Z 0 > + (N o KT l + N 1 KT o )(N 1 Z N + T 1 Z I ) 


1 3 


MN 


2 3x 


N l KT o n + Z o + N o Z N } + N o KT l n + z 0 + T oV 


- V 


3V 1 
1 3x" 


(5.9) 
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A - -It If N 1 KT 1 0 + Z „> + I N l KT o < T 1 Z T + N 1 Z N> 


+ f N 0 KT, (T^ + N,Z n ) + N, (T 1&r + N,^) 


- V 1 I? [ft 1 + Z o> « N l T o + W + ! V T o < N 1 Z N +T 1 Z T> 

+ + N 0 (N,e N + T 1 e T )_ 


r dV, 3V, 

- I < n 1 z n + T 1 Z T> n o KT o ST - I (1 + Z o>< N l CT o + W 5T 


av, 3 V-. 

N o ^ N l e N + t 1 £ T^ 3>T " N l e o 3x~ 


Next take 3/3t of (5.5) and eliminate N-j using (5.4) in both 
(5.5) and (5.6): 


/ 3V, \ 3T, 

v i + a k + *r- 


O 3x 


/ 3u l\ • 3V 1 

c € *-"o3rr dT i* s ir 


■ €A + e ( eN i + k ro 3 


Using (5.12) to eliminate T-j from (5.11) then yields the basic nonlinear 
fission-driven wave equation. 



3 2 V- 


3t 


, - V s 


3 2 V-i 
2 L 


3x^ 




mode coupling 



3A r 3£ 
3x + dt 



3 

3x 


( N l e + K R0 



(fission) (radiation) 


(5.13) 


where the sound speed V s is 

V = ! (g - c NJ 7 + a N | 

s l o' d o ) 

which readily reduces to equation (4.20) in section 4. 


(5.14) 


Now in lowest order (superscripts (o)) the general solution to 
(5.13) for the slab (which satisfies the boundary condition V-j = 0 at 
x = 0,L) is 


V 


(o) 

1 


z 

n=_oo 


i to t • i 
1 / e n sln k n x 



(5.15) 


(5.16) 


where 
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and 


,, _ un 

k n - IT 1 “n 


TrnV 

k n V s - TT 


* ■* * 
“ • W -n = W n ’ T -n = T n 


(5.17) 


iN 


W 


n V s *n 


(5.18) 


V s / aN c 
T n = 1 b 1 1 " UT 


(5.19) 


We must now derive equations for the complex amplitudes V through 
0(c). 


Multiple Time Scale Analysis 


Our objective is to find an asymptotic solution to (5.13) correct 
through 0 (c) and one that remains so through times of 0 (c" ). To do this 
we use a multiple time scale form of perturbation theory^ . Writing 
(5.13) as 


3 2 V 7 3 2 V 

— - - v* 

at 2 s 3x 2 


G L 


(5.20) 


we seek a solution of the form 


V-j = (t,ct,c 2 t, ..) + cvj 1 ^ (t,et, ..)+. 


(5.21) 
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where vj 0 ^ is chosen so that it reduces to (5.15) as C-> 0. Similar expan- 
sions are made for N-j and T-j and the boundary conditions must be satisfied 
through all orders. 

Two-Mode System 

The critical wavelength, 2-ir/k . t , above which instability occurs 
is fairly large so that it is expected that many systems will contain at 
most only a few unstable modes. In order to simplify the analysis we shall 
now treat the case in which only the fundamental, k-j = tt/L, is unstable and 
the second harmonic, !<£ = 2m/ L, is radiation-damped. Since higher harmonics 
(n >3) are progressi vely more heavily damped than the n = 2 mode we can 
neglect them and approximate the system by a 2-mode system with one mode 
unstable and one mode damped. 


( unstable 



h 

energy 

damped 


• 

^2 

^3 — 0 


1 

i 


Figure 5.3 


In this case 


/ \ 2 i 00 t 

= S V e sin k x 
1 n=-2 n " 


(5.22) 


and to 0(c) equation (5.20) becomes 
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where = !_ ( N j 0 ^ , vj 0 ^, T^). We must use the Ct dependence of 

to cancel those components of that would solve (3 2 /3t 2 - V|3 2 /3x 2 ) 

L resonant parts = These components, unless canceled, would give 
unwanted secular (t or x-proportional ) terms in V^. 


Equating the first term on the left of (5.23) to the resonant 
(o) 


terms of L v 'we have 


3U Hi) t . , 

0 . v n n sin k x = 
2i Z or - A - r - e n 


n 3 et 


Terms of of the type 


sin (k-j 2 x ) e 


11 W 1 ,2 t 


(5.24) 


where 


/be \ 3^ kM+M k JL 

\ d " a / 3x " d 3x 3t " d 3x 


3 Z T i 


J l e + K R0~ 


3x 2 J 


(5.25) 


Evaluation of Resonant Terms in Eg. (5.24) . 

By inspection of (5.25) and (5.8) - (5.10) for ip, <j>. A, we 
see that the following terms must be evaluated: 


3x : 


( N lV’ 3x3t ^iV’ 3x3t 


T 


2 

1 


3 

3t 



3N l\ 9 
9 xj’ 9t 
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(5.26) 

3_( v !M i_L Hi) 

3x \ V 1 3x / * 3x \ V 1 3x / 

3 _( t !M W n !M 

3x y 1 3x ) * 3x \ 1N 1 3x / 


where the small density derivatives Z^, will be neglected.* 

The resonant terms will be extracted from the first of these in detail and 
the result for the rest simply listed since they all involve similar cal- 
culations . 


Thus , 


p,2 2 i(w +w )t 

(N-.V-.) = I N 1/ e cos k x sin k x 

8x 2 1 1 3x 2 n,m=-2 nm n m 


Noting cos x sin y = [sin (x + y) - sin (x - y)]/2 this becomes. 


1 1 '( W n +u rJ t 

= - l X NV e n m [ (k + k ) 

2 n m L v n nr 


n = 2, m = -1 
contributions \n = - 1 , m = 2 

to sin k-jX \ n = -2, m = 1 

n = 1 , m = -2 


n nr 



Typically (see section 6) Te-j. ~ 12 * 1 0“ 1 1 , Ne^ ~ 2 • 10" 11 , 
TZ t ~ 2.6, NZ n ~ 0.26. 
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contribution ( n = 1, m = 1 
to sin k 2 x ( n = -1 , m = -1 


Ot 


n = 1 , m = -1 \ sin k 2 x e 


nonresonant 


Extracting the resonant terms leads to. 


i co *i t 


ax* 


1 I UJn U 

(N-jV-j) - j (^2^-1 + W -1^2^ e ^1 S1 ' n ^l x 


iw n t 


1 - ILU-j 

+ Y (N 2 ^-] + Wi t/ 2 ) k^ sin k-jX e 


1 / iwot -iw ? t \ 

e c - W_ 1 l/_ 1 e ) k| sin k 2 x 


-I lW,t k o i COo t 

j (^2^-1 + w _i ^ 2 ^ e k l sin k l x " 2^ e S1 ' n k 2 x H 


The remaining terms of (5.26) can be listed: 


1 t 


8x3t 


• I W-J l. 

( T-j ) -> - 7 >-w-]k.| (^ 2^-1 + W - 1 ^V e S1n k l x 


^ i oj?t 

2 e sin k 2 x + c.c. 


i w-j t 


i w 0 t 


3x3t ^1^ ~ ^1 ^ l ^ 2 ^ 1 ^ sin k 1 x — 2 T-jW 2 k 2 e sin k 2 x 


(5.27) 


c.c. 


(5.28) 


(5.29) 


+ c.c. 
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3T 


1 \ i (a)"| t 

i ( k i A / 9 T_ 1 - k 9 W_ l T 9 )ta- I e sin k-jX 


3 t V 1 3 x / 2 ^r 2 '-l 2 - i 2 ' 1 


• i cop t 

W-jT-j e k ^2 sin k 2 X + c.c. 


(5.30) 


a / 3N,\ _• i^t 

3t \ N 1 3x/ ' 2 k l u l W 2 W -l e sin k] x 


• 1 CD ? t 

j k-jOOgW^ e sin k 2 X + c.c. 


(5.31) 


3 / v 3V i 


Stnix 2 2-1 


i k -i co-j i ojn t i k -i cop i cop t 

l/ 0 l/ n e sin k-j x + — e sin k 2 x 


+ c.c. 


(5.32) 


3 /», 3N > 


3x VI 3X 


k? iWnt 

2~ ~ 2 V 1^2^ e s ^ n k ] x 


1 U)/p "t 

l /-j W-j k=j e sin k 2 x + c.c. 


(5.33) 


Uv, 3T ’ 


3x VI 3x 


k? i(o,t 

2 1 (t/ 2 r _] " 2t/ _i T 2^ e sin k l x 


ico 0 t 


Jrp l 

- IZ-jT-j k^ e sin k 2 x + c.c. 


(5.34) 
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3x 



k 2 

K 1 


(V-i 


i(i), t 

2T_ i ) e sin k-jX 


iout 

k| e sin k 2 X 


+ c.c. 


3 

3x 




k 2 

T~ ( N 2 U ~] 


id), t 

2W_^ 2 ) e sin k^x 


i d)pt 

- W, l/, k| e sin k2X + c.c. 


(5.35) 


(5.36) 


We also require for the last part of (5.13), 



+ K R0 



10), t 

eW l + K RO k l T l ^ e sin k l x 


2k i (~ eW 2 + 4K R0 k l T 2^ 


i 0 ) 2 1 

e sin kgX 


(5.37) 


The remaining algebra is tedious but straightforward. It con- 
sists of finding equations for SL^/Set and ai^/Sst from (5.24) using 
(5.25) for L, (5.8) - (5.10) for ip, A,cj>, together with (5.27) - (5.36) 
and the relations (5.18) and (5.19). The resulting equations become. 


31/ 1 

j Set = Y l l/ 1 + “V-l 

I 31/p 

3et r 2 V 2 + al/ l 


(5.38) 
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where 


Yi = 


2V. d 
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eN 0 V, 

2. . 1(2 K 

L V S K 1 N R0 b 


0 -? -w- 
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aN, 


'‘l 


k 2 - k 2 
k crit k l 


y 2 2V s d 
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V_ “ 4 k l K R0 ( ] w2 
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4k 2 
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k l j K 

2w ] d |2M 


e T Z T + (1 + Z ) (e - T e T \ 
ool o o oT) 


N„£n 


aN, 


KZ v 3 0 T 1 ' "° 


v s d 


aN, 


‘ T s Mb 2 


1 -ZT -i--T N o v s 1 ~ (Mb + t E r> < 5 - 40 > 


n 




N b 


- Oo + t(' + Z o> KT o]2? 


Nonlinear Behavior of Unstable Waves and their Limiting Amplitude 

We now have a simple pair of equations describing the nonlinear 
evolution of an unstable wave of velocity amplitude l/-j . Equations (5.38) 
contain the destabilizing effect of the fission power density together with 
the stabilizing effects of radiation diffusion and mode coupling to the 
"nearest" damped mode in wavenumber space. The smallness parameter c can 
now be set equal to unity and equations (5.38) for the wave amplitudes 
become 
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at 

ai / 2 
at - 


Y 1 l/ 1 + al/ 2 l/_-j 


+ al/f 


(5.41) 


(5.42) 


Now by adding I/, times (5.41) to t/ -j times (5.41) and noting 
from (5.17) that l' 1 = -^ (and carrying out a similar operation for 
(5.42)) these equations become 


ft hi i 2 * ’ 2y i I ‘'l 1 2 - 1/^*3 

(5.43) 

| ¥ ii/ 2 r--2Y 2 i(/ 2 r+ o[i/ 2 i/f + i/ 2 i/* 2 ] . 

i 01 p 

Writing the complex amplitudes as 2 = | l/, 2 |e 5 , these reduce 

further to 


|t l^l 2 = 2 y 1 |1/ 1 | 2 - 2a | ^-j | 2 1 ^2 i cos (0 2 - 26,) 

■gt I ^2 1 = " ^2^2^ + 2ot 1 1^-j | 2 1 1^ 2 1 cos (0 2 ■ 20, ) 


(5.44) 
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Suppose now a wave l/-j starts to grow out of the thermal noise 
and becomes of large amplitude. Since thermal fluctuations l ^ exist with 
all possible phases e 2 , the mode l/ 2 which begins to grow due to mode 
coupling will be that for which the mode-coupling term in the l/ 2 equation 
is largest, i.e., the mode 1 / 2 with a phase such that cos(0 2 - 20-j) = ±1. 
The limiting amplitudes at whi ch a_ steady state is reached fol low from 
(5.44) as 



^2 

a 2 

*1 


(5.45) 


where y-| , y 2 and a are given by (5.39) and (5.40), and the negative sign 
(cos(0 2 - 20-j) = -1) has been taken since a is negative. 
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Energy Density of a Nonlinear Wave 


From (2.31) the averaged wave energy density follows as 
W wave ■ <? MN o V l 2 + i N l KTtO + Z 0 ) + f (N, KT 0 + N 0 K Tl )T,Z T 


+ N-j T 1 e T > 


(5.46) 


i co-i "t 


"i 00 t "t 


1 UJ-j U “ I UJ-| L 

Now V-j (fundamental ) = l/-,e sin k-^x - e sin k-jX 


= 2 1 1/^ | sin k-jX cos (0 + o^t) 

* i i i 0 

where we used l/_-j = and set l/-j = |l/^|e . Thus, 


<vf> 

= l^l 2 



<N ] T 1 > 

- i^ 2 b° (’ - v ,° 


v s 

< T 1> 

-l^l 2 


where <> denotes a space and time average. Finally, using the relations 
(5.18) and (5.19) gives for the wave energy. 
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wave 






K 




( 5 . 47 ) 


It is useful to compare this with the thermal energy and define 

the ratio 

W 

r _ wave 

!vV 1+z o> 

2 1 U-j | 2 

3I W 1+ V 

( 5 . 48 ) 



-) ♦ I V Z T - 6 
v* / 2 0 T b s V 



wave energy 
thermal energy 


A similar expression exists for the second harmonic, etc. 


In the next section some numerical examples will be given based 
on the various formulas derived in this and the previous section. 



6. NUMERICAL CALCULATIONS USING PLASMA PARAMETERS CORRESPONDING TO THE 
NASA-LEWIS GAS-CORE DESIGN ENGINE 


In this section we evaluate the various expressions for the growth 
rates, limiting amplitudes, etc. for a uranium plasma with characteristics 
correspond!' ng to the NASA-Lewis gas-core design engine. The parameters were 
kindly supplied to us by R. Ragsdale of NASA-Lewis and are listed below: 



V. 


P fiss 


M 

P 

N. 


total U 


Tot 


50,000 °K 

1.3 * 10 10 ergs/cm 3 /sec = power density 
4.45 *10 6 cm 3 = U-plasrna volume 

2.7 • 10 4 gm. = total U-mass in reactor 

6.07 ‘ 10~ 3 gm/cm 3 = U-density 

1.5 • 10 19 cm" 3 = U number-density 

400 atmospheres =4 • 10 8 ergs/cm 3 = pressure. 


( 6 . 1 ) 


f 8) 

For the Rosseland opacity we use the calculated data of Parks et al v . 

The mean Rosseland opacity is 

k R = 2 * 10 4 cm 2 /gm. for P = 4 * 10 8 , T = 5 ■ 10 4 °K v 

i.e.j since p =6.07 • 10” 3 cm -3 it becomes 

k R = 1.2 * 10 2 cm' 1 (6.2) 

when expressed as an inverse length in the plasma. 


The corresponding radiation diffusion coefficient follows as 
K r = - 3.2 * 10 8 erg cm" 1 sec -1 deg" 1 


(6.3) 


where o is the black-body constant. 
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Evaluation of mean Charge, Z, mean Ionization Energy, e, and their 
Partial Derivatives e T , Zy, . 


The following quantities can be read off directly from the 
UPLAZ-2 graphs shown in Figs. 3.1 and 3.2 for a plasma at 50,000 °K and 
pressure of 400 Atm.: 







T 

Te 

TZ 


= 48 eV e 7.7 • 10 _11 ergs 
= 3.8 

= 2.4 * 10 -1 S ergs/deg. 

= 5.37 • 10" 5 deg -1 
= 1.2 * 10 _10 ergs 
= 2.68 


(6.4) 


The derivatives Z^ (= 3Z/3N) and e N are also of interest. They 
involve slightly more calculation since the graphs in Figs. 3.1 and 3.2 
give Z, e as functions of P, T instead of N, T. The conversion can be made 
by noting 


3_Z _ 3_ 

3N 3N 


z (nkt(i + z), t) 


3Z 

3P 


KT( 1 + Z) + NKT 


3Z 

3N 


from which it follows that 

3Z 


NZ N = 


N ~ (1 + Z) KT 


1 - NKT 


3Z 
3P J 


1 - 


3Z 

~3P 

P 3Z 
TT+I) 3P. 


N £., = P 


3e 


N ' 3P 


1 + 


3Z 

3P 


(’ + z - p i) 


(6.5) 
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To obtain e p arid Z p we make use of graphs of e(P) and Z(P) shown 
in Figs. 6.1 and 6.2 which were constructed from Figs. 3.1 and 3.2. These 
gave 


P e p = 12.6 eV = 2.02 • 10' 11 ergs 

P Z p = .25 

so that 

N Z N = .264 

N e^ = 2.13 • 10 _11 ergs. 


( 6 . 6 ) 


(6.7) 


We are now in a position to calculate various quantities of more direct 
physical interest. 

Sound Speed 

The sound speed as given by (4.20) is 


V s ■ 


5KT 

0 
3M 


, 3N e M 3 e t (1 + Z + N Z.,) 

| (1 + z 0 ) - ,^-N + T' o ON' 


5KT 


5K n + Z n + T o Z t) 


(6.8) 


K (in- Z 0 + T 0 Z t ) 


Noting KT Q = 7 * 10" 12 ergs, and 

KT P 

M p'TT +' Z) » (6.9) 
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and using the values (6.1) , (6.6) , (6.7) , we find 


V s - 2.73 • 10 s cm/sec. 


( 6 . 10 ) 


Critical Wavenumber 


The critical wavenumber for instability is given by (4.29). To 
evaluate we need the quantity. 



(a y 

v M h 


1 + z 


N o £ N 

KT 




'T 


1 KfT 


+ Z o + W 


It then follows from the above numerical values that 


k crit 


0 + Z o + W p fiss 

r 3 + , £ T . I 

) 

[2 K(1 + Z 0 + I 0 Z T )_ 

( T k ro 

", v N o £ N 3 N o Z N 

1 L KT 2 J 

1 


28 4 P ' " 2 
p fiss 


1 .35 T K 


R0 


= 1.32 * 10" 1 


( 6 . 11 ) 


( 6 . 12 ) 


The critical wavelength above which waves are unstable is thus, 


^crit 


2tt 

crit 


47.5 cm. 


(6.13) 
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Asymptotic Growth Rates 

For small k the growth rate in (4.28) becomes 

wj = 4.1 * 10 3 k 2 / 3 (6.14) 

At the other limit, xl.e., near marginal instability, we note 

Wj = 1*25 * 10 2 (k 2 rU - k 2 ) (6.15) 

where we used V 2 - A/M = 4.5 * 10 9 . 

An estimate of the maximum growth rate is 

“I Max 3 '- 2S • 102 k crit 3 2.18 sec -1 (6.16) 

These various values of Uj are plotted in Fig. 6.3. 



Figure 6.3. 


Growth rate of unstable waves as a 
function of wavenumber. 



We see that the growth rate is rather slow, i.e. 3 it would typically 
take a few seconds for the unstable waves to e-fold up in amplitude to their 
saturation value. This is to be compared with the real part of the frequency, 


which is of order (not exactly equal to) ^ .^.-V = 3.6 


10 4 sec -1 


Nonlinear Amplitude 

The limiting amplitude resulting from instability follows from 

(5.45) as 


(Vz )* 5 

■Rv: ' 


'RO 


(••?> , 


K >1 


L 


( k crit 


k l / \ k f 


k 2 ) 
crit /j 


(6.17) 


where 


{ ) = 


IM V7 L £ o T o Z T + ^ + Z > (e o - T o £ T> 


ILEr / 

K 7 Ml -2-1 h 

1 S Mb 2 


aN. 


*1 


V s d 


T K o V s 


J)( Mb + I £ t) 


N b 
o 


2V 


s L 


£ 0 + 1 n + z ) KT o 


We require the following quantities: 
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© » 

0 


© a 

- 5.83 • 

00 
t — i 

1 

o 

© a 

- 1 .98 - 

1 5 

10" 

©a 

- 2.77 * 

10' 16 

© a 

- 2.88 • 

10" 15 


Thus 


N V c 
o s 


{} s -5.14 • 10' 


■i s 


from which it follows that 



1.67 • 10“ 3 






h 


Ratio of wave energy to thermal Energy 


This ratio, which is given by (5.48) is 


9 NJW 2 
D _ 2_ r 2 _0 s 

K 1 3 1 P 


2 N Mb 


aN o\ 


V s 


7 


+ ? kz t 


MV 2 

s 

M 2 b 2 


1 - 


aN c 


where for the plasma parameters listed in this section the terms in the 


(6.19) 
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square bracket are 


[ ] s .5 + .418 + .003 = .921 . 


Thus 


I wave energy 
\ thermal energy 


0.675 



where |l/,|/V s is given by (6.19). 


( 6 . 20 ) 


Concluding Numerical Comments 

The above formulas (6.19) and (6.20) apply for a slab in which 

a single mode of wavenumber k-j is unstable due to coupling to a damped 

mode kg. Now recall that k-j = tt/L must satisfy k^ < k < k 2 = 2k-j 

in order for our two-mode analysis to be correct. Choosing for example 

kl = 2 k crit^ corresponding to a plasma slab of thickness 

L = 36 cm.) gives [ l/-j | /V s ~ 2 • 10 " 4 and R ~ 3 • 10~ 8 . 


The waves that result from the instability are thus seen to be of 
negligible amplitude in this case. This is because they have a relatively 
slow growth rate (wj/kV s ~ 10 4 ) and lose their energy via mode coupling to 
the radiation damped modes. For a thicker slab of plasma several modes 
become unstable and the algebra of the analysis becomes more complicated. 
However this is unlikely to change the above conclusion provided only a few 
modes are unstable, i.e. } provided L does not exceed a couple of meters. 

For example one would expect the energy in the fundamental to increase as 
L 5 / 2 (assuming a Kolmogoroff spectrum) which would increase the pressure 
fluctuation level to R ~ 10“ s for a slab of thickness 2.4 meters. 
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7. DEPENDENCE OF THE FISSION POWER DENSITY ON THE FLUID VELOCITY OF 
THE URANIUM, PLASMA AND ITS ROLE IN THE INSTABILITY OF SOUND WAVES 

The density of thermal neutrons within the reactor cavity consists 
principally of neutrons returned to the cavity from the surrounding moderator. 
Due to the open nozzle of the rocket chamber there will be an effective aver- 
age drift of the neutron gas out of the nozzle along with the hydrogen propel- 
lant. For purposes of calculation we shall assume that the thermal neutron 
flux in the rocket chamber can be represented by a Maxwellian distribution. 


A/ 

(2irv n ) 3 / 2 


(V - 


exp - 


V 2 




(7.1) 


where V n is the drift velocity towards the nozzle and v n the neutron thermal 
velocity. The space dependence of W or 1 will be neglected for the present, 
but should be discussed in a later study. 


Similarly the uranium plasma can undergo fluid motions with 
velocity V_ within the reactor cavity and the appropriate distribution function 
for the uranium nuclei is 


f 


u 


N 

(2tv u ) 3 / 2 


exp - 


(v - V) 2 



(7.2) 


Now due to the dependence of the neutron-induced fission cross- 
section on the relative velocity between a neutron and uranium nucleus, a 
coupling can occur between the plasma motion and the fission power density . 
The purpose of the analysis in this section is to examine the consequences 
of this coupling for the question of instability. 

The general expression for the fission power density can be written 
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P fiss = // d ^ d E Iv " v‘|Q f n (v) f u (v‘) a a (|v - v'|) 


(7.3) 


where Q is the kinetic energy in fission fragments in the plasma rest-frame 
(recall that in the energy equation we took the moment using i M(_v - _V) 2 
so that Pf.j ss is the rate of increase of thermal energy in the uranium plasma). 
The integral involves the Maxwellian distributions above in which it is of 
interest to note that the thermal velocities of the neutrons and uranium 
nuclei are comparable and the drift velocities and _V are subsonic. 


The integral required depends on the neutron absorption cross- 
section ct . This cross-section has a large number of resonances and it is 
usual to model these by a single resonance level using the Breit-Wigner 
formula (see Blatt and Weiskopf^, 



where 


E 

E 


o 

r 


and 



energy of resonance model for nucleus, 

kinetic energy associated with relative velocity 
between neutron and nucleus 

radiative capture and neutron widths 


(7.4) 


Now the energy E Q of the lowest neutron resonance is about 1 eV. 

In this analysis we will assume that the thermal neutron temperature is below 
10 4 °K in which case o & becomes approximately constant. This greatly simpli- 
fies the evaluation of (7.3) which then becomes. 
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P fiss (z-v*)-/* (a™;)*/* 

a 0 N N 
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Noting 


2v 


2v: 


/dvjv - A j e n = v< ^ v ^ y2 + ^ 2 ) e 0 


we find 


P fiss 


1 


WNa 0 - 

3 (2ir) 3 ^ 2 v 


rl 6 vj 


'-HI 


2 + 3v 


u 2 } 


Since v n , v u » j | , |Vj the expression for P^. can be written 


P_. = P-. + V * P' . 

fiss fiss,o — — fiss ,o 


where for the above assumed cross-section 


P‘ 


3P-. 
fiss ,o 


— fiss.o 3V 


2V P_. 

— n fiss,o 


l Vf° 3 (2v n 2 + v n 2 ) 


(7.5) 


(7.6) 


(7.7) 
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It should be emphasized that the general form (7.6) remains correct for the 
more complicated functional dependence of a as given by (7.4), but PL. 

Q I I o j j U 

would become then more complicated. It is interesting to note that a linear 
term in V_ only occurs if ^ 0, i.e. 3 the neutron population is drifting 
through the uranium population. The instability that results (as will be 
seen in the next section) is analogous to the two-stream instability of 
longitudinal waves that occurs in plasmas for which relative drift between 
species is present. 

Dispersion Relation for Sound Waves including the Velocity Dependence 
of the Fission Power Density. 

The fluid equations and linearization proceed as in section 4. 
Retaining the extra term (last term of (7.6)) the linear system of equations 
becomes , 


9N-j 

3t~ 


+ N o 



0 


(7.8) 


W, 1 3 

at + MN 0 3x 


N l KT o + Z o + W + N o KT l 0 + Z o + T oV 


= 0 (7.9) 


3 3_ 

2 at 


N l KT o n + Z o + N oV + N o KT l 0 + Z o + T oV 


at 


+ N 1 c q + N Q N 1 e N + 


+ 1 0 - z o> N o KT o Ir • ii 


(7.10) 


N 


+ w e — • V, 

‘Vo ax -i 


1 


N Q P fiss,o + -1 ' — fiss,o + K Ro v2T l 



As before we set 


(NpV-^T-,) - (N 1 ,V 1 ,T 1 )e l( ^" wt) , 

and eliminating N-j the transformed equations become 


(7.11) 


A BT 1 

+fek(k'V 1 ) + = 0 

0 


(7.12) 


- N ( i A + e + N \ + Pfis - i . v° + | (1 + l) N KT ft + N e J 
o \2 o o N I to 2 ooo o o J, 


k-V, 


^•1 ' — fiss.o ‘ *1 ) _ "V? u 1 Vt) ‘ 1 * “'Ro 


+ i V, • Pi, + T, -u iB + N 


( 


k 2 K 0 J = 0 
2 


(7.13) 


Eliminating T-j from (7.12), 


u -1 + Mw -1 ) — MN q 


( +i ^l *^fiss,0 


= 0 (7.14) 


Now the waves are longitudinal, i.e k_ is parallel or antiparallel to _V -j . 

Choose V_-j as x axis and set k = (k,0,0) (k has a sign). Then 

\ 

-1 * — ■ fiss o^— ’ -1 = P fiss o cos where 0 is the angle between _V-j 

and f.'f-j ss 0 - Taking k/ of (7.14) then leads to the dispersion 
relation 
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Vt ) • 1 
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RO 


k 2 B 
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■ N o (f* + e o + Vn ) 


+ l — D s . s . »° +1(1+ ZJN KT + N 
to 2 ' o' o o 


Mff P fiss,o cos 6 0 


This is again a cubic in to: 
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The sound speed .V can be defined as in section 4 and the cubic 


rewri tten 


Bk 


(f B + N o £ t) w ( w2 “ k 2 V s) + w2ik2|< Ro " iw P - cos 0 


MN q fiss,o 


Ai k 4 K 


Ro k 2 Bi 


M 
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i v iN Q fiss.o 



(7.15) 


(7.16) 
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Marginal Instability Condition 

For the case that the imaginary terms are small in (7.16) we can 

wri te 


01 = 0)0 + io) T = kV„ + 01 t 

R I si 


(7.17) 


and derive the marginal instability growth rate by linearizing in oij. 
follows from (7.16) that. 


It 


oi. 


MN ^ P fiss,o + V s P fiss,o cos " k2V s K Ro + 


Ak 2 K 


ro 


2N e T 

3 V H B + -F 


(7.18) 


The unstable wavenumber range now becomes. 




k I < k 


crit 


8 ( p fiss.o + Vfiss.o cos " 
ho m o k ' F 


(7.19) 


Small k Limit 

In the small -k limit the velocity-dependence of P^.. does 
not change the result found in section 4, so we again have 


(yj-t-i.) k 2/3 
2 K 


2 B P 


-i i/s 


fiss,o 


2N e T 

3M N o (B + -§ 


01 


(7.20) 



' General Comments 


The growth rate (7.18) near marginal instability behaves as if the 
fission power density was replaced by 


P,. + V 

fiss,o s 


p 1 

flSS ,0 


COS 0 


The velocity dependence of thus has a destabilizing influence if the 

second term is positive and is stabilizing if it is negative. 


Now from (7.7) we see that 


VP' „ cos 0 = 
s fnss,o 


2V s V n Pfiss.o cose 
3 < 2v n + v 7> 


(7. 


Thus the forward cone of waves (along V^) is unstable with growth rates 
that are largest for 0=0. However, since V s ~.v , 



and V n < v n , the absolute change in growth rate for the unstable waves that 
derives from this effect is fairly small. 

Effect of the fluid velocity dependence of the Fission power density 
on the Nonlinear behavior of Unstable Waves 


In this section we carry through the extra term, V_ • PL. 

I* "1 S S j 0 
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in the nonlinear multiple-time scale analysis of section 5. The nonlinear 
sound wave equation (5.13) then becomes 






mode coupling 



3A 

3x 


+ £ 


3«J> 

at 


_ b 3 
e d 3x 


N l e+ V 1 P fiss,o + K Ro 


S 2 Tn 


ax' 


(fission) (fission (radiation) 
Doppler 
effect) 


(7.22) 


If we make the same expansion as before, with a zeroth-order expression 
for V.j , 


co it 

vj 0) = Tj ^ e sin k n x ’ (7.23) 

n=-oo 

i co -j t 

it turns out that terms of the form sin k-jX e occur in 0(e) throughout 
(7.22), except for the extra term involving Pi- . This extra term only 

T I bb ,U _• i 

1 COt t 

contains first harmonic components of the form cos k,x e . Multiplying 

/ L * 

dx, the effect of the 
o 

P|. term therefore averages to zero, i.e. 3 for a bounded geometry it has 

T 1 SS )0 — 

no effect on the nonlinear growth of waves. 


The reason for this can best be seen by including the Pi, 

I I Oj j U 

term in the zeroth-order equation and considering the equation 
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3 2 V, 


3t‘ 


9 2 V n 
/ 2 L 


+ a 


3x‘ 


3V n 

dx 


= 0 , 


(7.24) 


where a = b P' f . /d, for the slab geometry. Substituting traveling wave 

1 ! b o j U 

solutions gi^-kx) gives a dispersion relation 


3 * 


w = ± (k 2 V 2 - iak) 2 = ± kV 1 


iak 


2k 2 V 2 


(7.25) 


so that a solution of wavenumber k satisfying (7.24) is 


at 


at 


V, S A e 2V e i(kVt - kx) + B e" 2V e - i(kVt+kx) 


(7.26) 


It is clear from (7.26) that the right traveling wave is growing 
in amplitude (a > 0), and the left-traveling wave decaying. Now for a 
bounded geometry with perfectly reflecting boundaries a standing wave con- 
sists of a superposition of a right and left traveling wave. As the wave 
traverses the slab thickness L to the right it grows in amplitude. Then, 
after reflection, it decays an equal amount on traveling to the left. The 
average effect is a zero growth for a standing wave in a slab. 

In the actual case of a gas core nuclear rocket engine, the only 
waves which can propagate freely out of the U-plasma without reflection, 
are those propagating out of the nozzle. These will undergo at most only 
a few e-foldings in the course of their transit across the U-plasma and are 
not expected to be of importance. 
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8. COMMENTS ON THERMAL FLUCTUATIONS IN A NEUTRON GAS AND IN THE DISTRIBUTION 
OF FISSION EVENTS IN A U-PLASMA 


The main point of the study described in this report is to calculate 
the behavior of those unstable waves in the system that may become of large 
amplitude. However it is worth noting that there is of course an irreducible 
thermal fluctuation level in such a plasma. These thermal fluctuations pro- 
vide the initial amplitude from which unstable waves may grow. Fluctuations 
in the fission power density, although small, will be larger than is normally 
the case due to the large amount of energy (200 MeV) released per fission 
event. 

Consider first the fluctuations in the noninteracting neutron gas 
of N neutrons in a box of volume V. The neutron distribution function can be 
written 


F = ~ . Z " X i ( t)l ’ V 

o i=l 

( 8 . 1 ) 

= Ai 0 + = V io = constant. 


'f* h 

where X_. , are the position and velocity of the i U1 neutron. Various 
statistical quantities can be calculated using an ensemble function for 
the initial particle parameters normalized so that 





dX N0 dV N0 


1 


and the s-particle distribution is defined by 


( 8 . 2 ) 


N,0 


f s ^10’ •*’ V s0 ,1: ) 


yS f °N dX s + l,0 ••• dV i 


(8.3) 



The correlation function of interest to us is that corresponding 
to fluctuations in the neutron number density, i.e. 3 

/dv, dv 2 <F(l,t 1 ) F ( 2 , t 2 ) >' d lz dX 10 •** dV N0 FF 


N 2 V 
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dV io f l ^ V io^ 6 ^-l “ ^2 " -io ^1 “ 


+ "fe 1 ) /dV • dV- f, (V. ) f , (V. ) 
W 2 V 2 J -TO -jo '1 —io' 1 v -jo' 


i + _1 1 

W- / j. j. \ 3 ' 1 l ti - t 


x-i - x Q \ <6N <5W > 

f, i =!— ±1 = i + — - — — 


o (t 1 -t 2 ) : 


'1 u 2j 


W o 


(8.4) 


where W Q is the average neutron number density and 6N Q the fluctuation in 
number density. 

Assuming a thermal Maxwellian distribution for f-j (normalized 

to unity) , 


gives 



1 

(2 irv’) 1 / 2 


exp 



(8.5) 
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(x-| - x 2 ) 2 


(2ir) 3 /2[v 0 (t r t 2 )] 3 


<6W q (x 1 ,t 1 )«SW 0 (x 2 ,t 2 ) > 


exp 


2v o (t rt 2 )“ 


(8.6) 
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where v is the neutron thermal velocity and N their average density, 
o o 

For fluctuations at a single point (x_-| = we have, 


6W 0 (l,t 2 )> 


1 

W q (2ti) 3 ^ 2 [v 0 (t r t 2 )]> 


(8.7) 


We see that this is small except for very fine scale fluctuations with 
periods of order 1/W Q 1 / 3 v Q = time it takes a neutron to traverse a mean 
spacing between neutrons. 

Fluctuations in the Number Density of Fission Events and in the Fission 
Power Density 

An approximate estimate of this can be obtained as follows. 


Write 

p . 

N^ss = no - of fissions/cm 3 /sec = ■■ f q -— . (8.8) 

Thus in a volume 5, 3 and time t the average number of fissions is & 3 t . 

Assuming these fission events are uncorrelated, the fluctuation in this 
number is U 3 x N fl - ss ) 3s . The corresponding fluctuation in the fission power 
density can be written 


6P 

~P 


fiss 

fiss 


1 

P fiss 

L Q 



(8.9) 
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Now Q - 200 MeV = 3. TO -4 ergs, and if for example we choose 
^fiss ~ 4.10 8 erg/sec/ cm 3 , it follows that 

6P fiss „ 1 . (8.10) 

P fiss [1.6 • 10 12 £ 3 t f 2 

The amplitude of 6P^. for some length and time scales are shown below: 


a = 

10“ 2 cm. 

II 

O 

t 

N> 

sec. «P fiss / p f jss 

a 10" 

a = 

10“ 3 cm, 

T = 10 _1 

sec, 6 p fiss /P fiss 

a 10“ 


We see that the fission power density has a "grainy" fluctuation level that 
is appreciable for length scales of ~ 10“ 2 - 10“ 3 cm, and for the time 
scales above. For lower power densities this irreducible thermal fluctua- 
tion becomes relatively more important. However, since they involve short 
length scales it is expected that radiation diffusion will considerably 
decrease the corresponding pressure fluctuations in the plasma. 



9. CONCLUDING REMARKS AND RECOMMENDATIONS 


The overall conclusions we reached as a result of this analysis 

can be stated as follows: 

(i) A range of unstable sound waves exists in a fissioning uranium 
plasma in the wavenumber range k < k ^ (see 4.29). For para- 
meters corresponding to the NASA-Lewis design engine the critical 
wavelength above which waves are unstable is X .. = 2-rr/k . + = 

47.5 cm. Typical e-folding times are .5 sec. The mechanism giving 
rise to the instability can be described as follows: Consider a 

standing sound wave in a bounded region of fissioning plasma with 

a constant background thermal neutron density. In the wave com- 
pressions the fission power density, P^.. , increases due to the 

increased uranium density, and in the rarefactions P fl - ss decreases. 
This results in an increased pressure gradient associated with the 
wave which in turn leads to a transfer of fission power to the 
wave. It occurs because the wave compression tends to expand more 
rapidly than it was compressed. However, competing with this is 
the fact that radiation tends to transport the extra thermal energy 
out of the wave compressions. Radiation diffusion smooths out the 
temperature fluctuations of waves more rapidly the shorter their 
wavelength. This results in a critical wavelength, X below 
which waves are stable and above which they are unstable. 

(ii) Waves that are unstable grow in amplitude until limited by nonlinear 
effects. The mechanism that finally limits their amplitude is the 
mode coupling of unstable waves in the wavelength range X > X 

to radiation damped waves in the range X < X ... Our conclusion 
is that this mechanism "turns-off" the instability after the wave 



amplitudes have reached relatively small values for most systems 
of interest. Relative pressure fluctuations of order 
5P/P ~ 10" 7 to 10 -5 are calculated to occur for a slab of 
uranium plasma with parameters corresponding to the NASA-Lewis 
design engine. Consequently this instability is not expected 
to present a control problem for this engine. It could become 
a problem, however, if very large devices were designed. We 
estimate that the pressure fluctuation amplitude would increase 
as the 5/2 power of the radius of the uranium plasma core. 

o 

(tit) The dependence of the fission power density on the fluid velocity 
of the uranium plasma in a cavity in which the neutron gas has a 
drift motion out of the nozzle can contribute to instability. 

The mechanism is similar to that for the above instability and 
is analysed in section 7. It also gives rise to small -amplitude 
fluctuations in the reactor cavity. 

The role of spatial gradients in the neutron density and plasma 
density has not been included in the above analysis. Since these gradients 
may also give rise to instability we recommend that these be investigated as 
a next step in this analysis. Further, the interaction between different 
flow states of the uranium plasma in the cavity and the fission power should 
be analyzed ( e.g. 3 the hydrogen flow may cause the U plasma core to circulate 
in convective cells). A computer simulation of these reacting flow problems 
would be of value. 
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